#ifndef CPP_DSLASH_MATVEC_64BIT_H
#define CPP_DSLASH_MATVEC_64BIT_H


#include <cpp_dslash_types.h>
#warning "Using C stuff for 64 bit"


namespace CPlusPlusWilsonDslash {
  using namespace Dslash64BitTypes;

    // Gauge Matrix is 3x3x2 (Natural ordering)
    // HalfSpinor is 2x3x2 (spin, color, reim)
    inline 
    void su3_mult(HalfSpinor res, const GaugeMatrix u, const HalfSpinor src) 
    {

      int re=0;
      int im=1;
	
      // spin 0
      // col = 0
      // row = 0
      res[0][0][0] = u[0][0][0]*src[0][0][0]
	-u[0][0][1]*src[0][0][1];
	  
      res[0][0][1] = u[0][0][0]*src[0][0][1]
	+u[0][0][1]*src[0][0][0];

      // row = 1
      res[0][1][0] = u[0][1][0]*src[0][0][0]
	-u[0][1][1]*src[0][0][1];
	  
      res[0][1][1] = u[0][1][0]*src[0][0][1]
	    +u[0][1][1]*src[0][0][0];

      // row = 2
      res[0][2][0] = u[0][2][0]*src[0][0][0]
	    -u[0][2][1]*src[0][0][1];
	  
      res[0][2][1] = u[0][2][0]*src[0][0][1]
	+u[0][2][1]*src[0][0][0];

      
      // col = 1
      // row = 0
      res[0][0][0] += u[1][0][0]*src[0][1][0]
	    -u[1][0][1]*src[0][1][1];
	  
      res[0][0][1] += u[1][0][0]*src[0][1][1]
	    +u[1][0][1]*src[0][1][0];

      // row = 1
      res[0][1][0] += u[1][1][0]*src[0][1][0]
	    -u[1][1][1]*src[0][1][1];
	  
      res[0][1][1] += u[1][1][0]*src[0][1][1]
	    +u[1][1][1]*src[0][1][0];

      // row = 2
      res[0][2][0] += u[1][2][0]*src[0][1][0]
	-u[1][2][1]*src[0][1][1];

      res[0][2][1] += u[1][2][0]*src[0][1][1]
	+u[1][2][1]*src[0][1][0];
	  
	  
      // col = 2
      // row = 0
      res[0][0][0] += u[2][0][0]*src[0][2][0]
	-u[2][0][1]*src[0][2][1];
	  
      res[0][0][1] += u[2][0][0]*src[0][2][1]
	    +u[2][0][1]*src[0][2][0];

      // row = 1
      res[0][1][0] += u[2][1][0]*src[0][2][0]
	-u[2][1][1]*src[0][2][1];
      
      res[0][1][1] += u[2][1][0]*src[0][2][1]
	+u[2][1][1]*src[0][2][0];

      // row = 2
      res[0][2][0] += u[2][2][0]*src[0][2][0]
	-u[2][2][1]*src[0][2][1];
      
      res[0][2][1] += u[2][2][0]*src[0][2][1]
	+u[2][2][1]*src[0][2][0];


      // spin = 1
      // col = 0
      // row = 0 
      res[1][0][0] = u[0][0][0]*src[1][0][0]
	-u[0][0][1]*src[1][0][1];
	  
      res[1][0][1] = u[0][0][0]*src[1][0][1]
	+u[0][0][1]*src[1][0][0];

      // row = 1
      res[1][1][0] = u[0][1][0]*src[1][0][0]
	-u[0][1][1]*src[1][0][1];
	  
      res[1][1][1] = u[0][1][0]*src[1][0][1]
	    +u[0][1][1]*src[1][0][0];

      // row = 2
      res[1][2][0] = u[0][2][0]*src[1][0][0]
	    -u[0][2][1]*src[1][0][1];
	  
      res[1][2][1] = u[0][2][0]*src[1][0][1]
	+u[0][2][1]*src[1][0][0];

      // col = 1
      // row = 0
      res[1][0][0] += u[1][0][0]*src[1][1][0]
	    -u[1][0][1]*src[1][1][1];
	  
      res[1][0][1] += u[1][0][0]*src[1][1][1]
	    +u[1][0][1]*src[1][1][0];
	
      // row = 1
      res[1][1][0] += u[1][1][0]*src[1][1][0]
	    -u[1][1][1]*src[1][1][1];
	  
      res[1][1][1] += u[1][1][0]*src[1][1][1]
	    +u[1][1][1]*src[1][1][0];

      // row = 2
      res[1][2][0] += u[1][2][0]*src[1][1][0]
	-u[1][2][1]*src[1][1][1];

      res[1][2][1] += u[1][2][0]*src[1][1][1]
	+u[1][2][1]*src[1][1][0];
	  
	  
      // col = 2
      // row = 0
      res[1][0][0] += u[2][0][0]*src[1][2][0]
	-u[2][0][1]*src[1][2][1];
	  
      res[1][0][1] += u[2][0][0]*src[1][2][1]
	    +u[2][0][1]*src[1][2][0];

      // row = 1
      res[1][1][0] += u[2][1][0]*src[1][2][0]
	-u[2][1][1]*src[1][2][1];
      
      res[1][1][1] += u[2][1][0]*src[1][2][1]
	+u[2][1][1]*src[1][2][0];

      // row = 2
      res[1][2][0] += u[2][2][0]*src[1][2][0]
	-u[2][2][1]*src[1][2][1];
      
      res[1][2][1] += u[2][2][0]*src[1][2][1]
	+u[2][2][1]*src[1][2][0];



    }

    inline 
      void su3_adj_mult(HalfSpinor res, GaugeMatrix u, HalfSpinor src)
    {

      // spin=0
      // row=0
      // col=0
      res[0][0][0] = u[0][0][0]*src[0][0][0]
	+  u[0][0][1]*src[0][0][1];
      res[0][0][1] = u[0][0][0]*src[0][0][1]
	- u[0][0][1]*src[0][0][0];

      // col=1
      res[0][0][0] += u[0][1][0]*src[0][1][0]
	+  u[0][1][1]*src[0][1][1];
      res[0][0][1] += u[0][1][0]*src[0][1][1]
	- u[0][1][1]*src[0][1][0];

      // col=2
      res[0][0][0] += u[0][2][0]*src[0][2][0]
	+  u[0][2][1]*src[0][2][1];
      res[0][0][1] += u[0][2][0]*src[0][2][1]
	- u[0][2][1]*src[0][2][0];

      // row 1 
      // col=0
      res[0][1][0] = u[1][0][0]*src[0][0][0]
	+  u[1][0][1]*src[0][0][1];
      res[0][1][1] = u[1][0][0]*src[0][0][1]
	- u[1][0][1]*src[0][0][0];
      
      // col=1
      res[0][1][0] += u[1][1][0]*src[0][1][0]
	+  u[1][1][1]*src[0][1][1];
      res[0][1][1] += u[1][1][0]*src[0][1][1]
	- u[1][1][1]*src[0][1][0];

      // col=2
      res[0][1][0] += u[1][2][0]*src[0][2][0]
	+  u[1][2][1]*src[0][2][1];
      res[0][1][1] += u[1][2][0]*src[0][2][1]
	- u[1][2][1]*src[0][2][0];


      // row 2
      // col=0
      res[0][2][0] = u[2][0][0]*src[0][0][0]
	+  u[2][0][1]*src[0][0][1];
      res[0][2][1] = u[2][0][0]*src[0][0][1]
	- u[2][0][1]*src[0][0][0];
      
      // col=1      
      res[0][2][0] += u[2][1][0]*src[0][1][0]
	+  u[2][1][1]*src[0][1][1];
      res[0][2][1] += u[2][1][0]*src[0][1][1]
	- u[2][1][1]*src[0][1][0];
      
      // col=2      
      res[0][2][0] += u[2][2][0]*src[0][2][0]
	+  u[2][2][1]*src[0][2][1];
      res[0][2][1] += u[2][2][0]*src[0][2][1]
	- u[2][2][1]*src[0][2][0];
      
      
      // spin=1
      // row=0
      // col=0
      res[1][0][0] = u[0][0][0]*src[1][0][0]
	+  u[0][0][1]*src[1][0][1];
      res[1][0][1] = u[0][0][0]*src[1][0][1]
	- u[0][0][1]*src[1][0][0];
      
      // col=1
      res[1][0][0] += u[0][1][0]*src[1][1][0]
	+  u[0][1][1]*src[1][1][1];
      res[1][0][1] += u[0][1][0]*src[1][1][1]
	- u[0][1][1]*src[1][1][0];

      // col=2
      res[1][0][0] += u[0][2][0]*src[1][2][0]
	+  u[0][2][1]*src[1][2][1];
      res[1][0][1] += u[0][2][0]*src[1][2][1]
	- u[0][2][1]*src[1][2][0];
            
      // row=1
      // col=0
      res[1][1][0] = u[1][0][0]*src[1][0][0]
	+  u[1][0][1]*src[1][0][1];
      res[1][1][1] = u[1][0][0]*src[1][0][1]
	- u[1][0][1]*src[1][0][0];

      // col=1
      res[1][1][0] += u[1][1][0]*src[1][1][0]
	+  u[1][1][1]*src[1][1][1];
      res[1][1][1] += u[1][1][0]*src[1][1][1]
	- u[1][1][1]*src[1][1][0];

      // col=2
      res[1][1][0] += u[1][2][0]*src[1][2][0]
	+  u[1][2][1]*src[1][2][1];
      res[1][1][1] += u[1][2][0]*src[1][2][1]
	- u[1][2][1]*src[1][2][0];

      // row=2
      // col=0
      res[1][2][0] = u[2][0][0]*src[1][0][0]
	+  u[2][0][1]*src[1][0][1];

      res[1][2][1] = u[2][0][0]*src[1][0][1]
	- u[2][0][1]*src[1][0][0];
      
      // col=1      
      res[1][2][0] += u[2][1][0]*src[1][1][0]
	+  u[2][1][1]*src[1][1][1];
      res[1][2][1] += u[2][1][0]*src[1][1][1]
	- u[2][1][1]*src[1][1][0];
      
      // col=2
      res[1][2][0] += u[2][2][0]*src[1][2][0]
	+  u[2][2][1]*src[1][2][1];
      res[1][2][1] += u[2][2][0]*src[1][2][1]
	- u[2][2][1]*src[1][2][0];



    }
}

#endif

    



